tRNA methylation in fetal tissues

1 Intro

1.1 Aim

• Aquire tRNA gene methylation data for fetal tissues, in order to establish if tRNAs are hypomethylated in fetal tissues

(Yang et al. 2016) created “epiTOC” (Epigenetic Timer Of Cancer) and used DNAm data from fetal tissues:

Age-hypermethylated CpGs were filtered further, selecting only those with absent or low (beta <0.2) methylation across 52 fetal tissue samples encompassing 11 tissue types (cord blood (GSE72867), stomach, heart, tongue, kidney, liver, brain, thymus, spleen, lung, adrenal gland (Nazor et al. 2012)).

The methylation data takes the form of mixed 450k and 27k illumina bead-chip methylation array beta values.

(Nazor et al. 2012) Tissue specific methylation data is under: GSE30654, which contains samples other than the fetal samples so accessions for the tissues listed in (Yang et al. 2016) were selected, this totaled 51 samples from 10 tissue types.

The epiTOC paper referenced GSE72867 for cord blood this contained 450k array data for CD34+ cells from cord blood including 5 controls each with 3 technical replicates. As there was no additional information on which samples were used from this set all the control samples were used, taking the mean of the technical replicates.

Specific accessions for the selected samples were compiled into a list (see below)

2 Setup

2.1 libs

suppressPackageStartupMessages({
    library(tidyverse)
    library(broom)
    library(rvest)
    library(plotly)
    
    library(heatmaply)
    library(ComplexHeatmap)
})

nest <- nest_legacy
unnest <- unnest_legacy

if(!dir.exists("graphics")){
    dir.create("graphics", showWarnings = FALSE, recursive = TRUE)
}

2.2 data read-in

2.2.1 accession from samples used in epiTOC

epiTOCfetalTissues <- read_csv("data/epiTOC_fetal_DNAm_datasets.csv")

Parsed with column specification:
cols(
  stomach = col_character(),
  heart = col_character(),
  tongue = col_character(),
  kidney = col_character(),
  liver = col_character(),
  brain = col_character(),
  thymus = col_character(),
  spleen = col_character(),
  lung = col_character(),
  `adrenal gland` = col_character(),
  cord_blood = col_character(),
  cordBlood_groups = col_character()
)

epiTOCfetalTissuesL <- epiTOCfetalTissues %>% 
    gather(key = "tissue", value = "accession", stomach:cord_blood) %>%
    na.omit() %>%
    mutate(
        cordBlood_groups = if_else(
            grepl(x = tissue, pattern = "cord_blood"), cordBlood_groups, "NA"
        )
    )

2.2.2 tRNA data

tRNAprobes <- read_delim(
    delim = "\t",
    file = "../tRNA-GtRNAdb/450k_coresponding_hg19tRNAs.bed",
    col_names = c(
        "pChr","pStart","pEnd","probeID",
        as.character(
            read_delim(
                delim = "\t",
                file = "../tRNA-GtRNAdb/std_tRNA_header.txt",
                col_names = FALSE
            )
        )
    )
)

Parsed with column specification:
cols(
  X1 = col_character(),
  X2 = col_character(),
  X3 = col_character(),
  X4 = col_character(),
  X5 = col_character(),
  X6 = col_character(),
  X7 = col_character(),
  X8 = col_character(),
  X9 = col_character(),
  X10 = col_character(),
  X11 = col_character(),
  X12 = col_character()
)
Parsed with column specification:
cols(
  pChr = col_character(),
  pStart = col_double(),
  pEnd = col_double(),
  probeID = col_character(),
  tChr = col_character(),
  tStart = col_double(),
  tEnd = col_double(),
  tRNAname = col_character(),
  score = col_double(),
  strand = col_character(),
  thickStart = col_double(),
  thickEnd = col_double(),
  RGB = col_double(),
  blockCount = col_double(),
  blockSizes = col_number(),
  blockStarts = col_character()
)

tRNAge <- read_delim(
    delim = "\t",
    file = "../tRNA-GtRNAdb/GenWideSigWintRNAs.txt",
    col_names = "tRNAname",
    col_types = "c"
)

gwsBB6 <- read_tsv(
    "../tRNA-GtRNAdb/BB_GWS_tRNA.txt",
    col_names = "tRNAname", col_types = "c"
)

swsBB23 <- read_tsv(
    "../tRNA-GtRNAdb/swsBB23.tsv",
    col_names = "tRNAname", col_types = "c"
)

bltRNAs <- read_tsv(
    "../tRNA-GtRNAdb/tRNAs-in-hg19-blacklist-v2.txt",
    col_names = "tRNAname", col_types = "c"
)

2.2.3 Codon Usage data / related tRNA data

tRNA <- read.table("../tRNA-GtRNAdb/hg19-tRNAs-SeqStrPseu.bed") 

colnames(tRNA) <- 
c(
    as.character(
        read.table(
            "../tRNA-GtRNAdb/std_tRNA_header.txt",
            stringsAsFactors = FALSE
        )
    ),
    "pseudo","str","seq"
)

tRNA <- tRNA %>%
    extract(
        tRNAname,
        c("nmt", "aa", "codon"),
        "(\\w*)-?tRNA-(i?\\w{3})(?:\\w+)?-(\\w+)-",
        remove = FALSE
    )

palette <- c("#e6194b","#3cb44b","#424242","#0082c8","#f58231","#911eb4","#46f0f0","#f032e6","#d2f53c","#fabebe",
             "#008080","#e6beff","#aa6e28","#fffac8","#800000","#aaffc3","#808000","#ffd8b1","#000080","#808080",
             "#000000","#00FF00","#424242","#0000FF")

CodonUsageFreqHuman <- read.table(
    "../tRNA-GtRNAdb/CodonUsageFreqHuman.tab",
    header = TRUE,
    stringsAsFactors = FALSE
)

2.2.4 Methylation data

2.2.4.1 Web-scrape

methylationArrayData <- lapply(
(epiTOCfetalTissuesL %>% 
    pull(accession)),
    function(x){
        paste0( "https://www.ncbi.nlm.nih.gov/geo/tools/geometa.cgi?acc=",
                x,
                "&scope=full&mode=miniml"
        ) %>%
        read_html() %>% #html_structure()
        html_nodes("internal-data") %>% #pre,x = "pre"
        html_text() %>%
        read_delim(delim = "\t", col_names = FALSE) %>% #c("probeID","beta","detection-p")
        mutate(accession=x)
    }
)
methylationArrayDataS <- lapply(methylationArrayData, function(x) {x %>% select(X1,X2,accession)})
methylationArrayDataDF <- do.call(rbind, methylationArrayDataS)
methylationArrayDataDF <- methylationArrayDataDF %>% rename(probeID = X1, beta = X2)

saveRDS(
    methylationArrayDataDF,
    file = "data/epiTOCFetal_methylationArrayDataDF.Rds"
)

2.2.4.2 Load Local data copy

methylationArrayDataDF <- readRDS(
    file = "data/epiTOCFetal_methylationArrayDataDF.Rds"
)

3 Data Overview

3.1 N probes per sample

methylationArrayDataDF %>%
    group_by(accession) %>%
    summarise(n = n())

3.2 Beta density plots

ggplotly(
    methylationArrayDataDF %>%
        ggplot(aes(beta, colour = accession)) + 
            geom_density() + 
            guides(colour=FALSE)
)

4 Data Processing

4.1 Mean betas from technical replicates

NB keeps an aribtary accession number from each sample for cross referencing with tissue type data

meanBloodBetas <- epiTOCfetalTissuesL %>% 
    filter(!cordBlood_groups == "NA") %>% 
    extract(
        col = cordBlood_groups,
        into = c("sample", "replicate"),
        regex = "(Ctrl\\d+)_(\\w)",
        remove = FALSE
    ) %>%
    group_by(sample) %>%
    do(meanBetas = {
        data <- .
        methylationArrayDataDF %>% 
        filter(
            accession %in% (data %>% pull(accession) )
        ) %>% 
        group_by(probeID) %>% 
        summarise(beta = mean(beta), accession = data$accession[1])
    })

saveRDS(
    meanBloodBetas,
    "data/epiTOCFetal_meanBloodBetas.Rds"
)

meanBloodBetas <- readRDS(
    "data/epiTOCFetal_meanBloodBetas.Rds"
)

Unnest and replace data with technical replicates with the mean values

methylationArrayDataReptMeansDF <-
bind_rows(
    methylationArrayDataDF %>%
        filter(
            accession %in% (
                epiTOCfetalTissuesL %>% 
                filter(cordBlood_groups == "NA") %>%
                pull(accession) 
            )
        ),
    meanBloodBetas %>% 
        unnest() %>%
        select(probeID, beta, accession)
)

saveRDS(
    methylationArrayDataReptMeansDF,
    file = "data/epiTOCFetal_methylationArrayDataReptMeansDF.Rds"
)

4.1.0.1 Load Local data copy

methylationArrayDataReptMeansDF <- readRDS(
    file = "data/epiTOCFetal_methylationArrayDataReptMeansDF.Rds"
) 

4.2 Get tissue categories

methylationArrayDataReptMeansDF <- left_join(
    methylationArrayDataReptMeansDF,
    epiTOCfetalTissuesL %>% select(tissue, accession),
    by = c("accession", "tissue")
)

5 Get meth data from tRNA probes

tRNAprobesMethData <- 
left_join(
    tRNAprobes,
    methylationArrayDataReptMeansDF,
    by = "probeID"
) %>% 
    extract(tRNAname,
            c("nmt","aa","codon"),
            "(\\w*)-?tRNA-(i?\\w{3})(?:\\w+)?-(\\w+)-",
            remove = FALSE
    ) %>% 
    mutate(tRNAge = tRNAname %in% (tRNAge %>% pull(tRNAname)))

tRNAmeanMethByTissue <-
tRNAprobesMethData %>%
    dplyr::select(probeID, tChr, pStart, tStart, tRNAname, nmt, aa, strand, codon, beta, accession, tissue, tRNAge) %>%
    dplyr::group_by(tChr, tStart, tRNAname, nmt, aa, strand, codon, tissue, accession, tRNAge) %>%
    summarise(beta = mean(beta))

tRNAmeanMethByTissue %>% nrow()

[1] 5344

tRNAmeanMethByTissue <- 
tRNAmeanMethByTissue %>%
    filter(!tRNAname %in% bltRNAs)

tRNAmeanMethByTissue %>% nrow()

[1] 5344

5.1 Number of tRNAs covered:

NtRNAgenes <- tRNAmeanMethByTissue %>%
    ungroup() %>%
    select(tRNAname) %>%
    distinct() %>%
    nrow()
NtRNAgenes

[1] 150

5.2 Number of tRNAge genes covered

NtRNAgeGenes <- tRNAmeanMethByTissue %>%
    ungroup() %>%
    filter(tRNAge == TRUE) %>%
    select(tRNAname) %>%
    distinct() %>%
    nrow()
NtRNAgeGenes

[1] 17

meanMethBytRNAByTissue <- tRNAprobesMethData %>%
    select(probeID,tChr,pStart,tStart,tRNAname,nmt,aa,strand,codon,beta,accession,tissue) %>%
    group_by(tChr,tStart,tRNAname,nmt,aa,strand,codon,tissue) %>%
    summarise(beta = mean(beta)) %>%
    ungroup()

# tRNAs with mean beta < 0.2 in 1 or more tissues
meanMethBytRNAByTissue %>%
    filter(beta < 0.2) %>% 
    select(tRNAname) %>% 
    distinct() %>%
    nrow()

[1] 130

# tRNAs with mean beta >= 0.2 in 1 or more tissues
meanMethBytRNAByTissue %>%
    filter(beta >= 0.2) %>% 
    select(tRNAname) %>% 
    distinct() %>%
    nrow()

[1] 48

5.3 tRNA mean metylation by tissue plots

5.3.1 Heatmaps

meanMethBytRNAByTissueMat <- 
tRNAmeanMethByTissue %>%
    ungroup() %>%
    select(tRNAname,tissue,beta) %>%
    group_by(tRNAname,tissue) %>%
    summarise(mean = mean(beta)) %>%
    drop_na(tissue) %>%
    spread(tissue,mean) %>%
    na.omit() %>%
    column_to_rownames("tRNAname") %>%
    data.matrix()

meanMethBytRNAByTissueMat %>% nrow()

[1] 115

col_fun <- circlize::colorRamp2(
    c(0, 0.5, 1),
    c("#ffff00", "#19ff00", "#0033cc")
)

meanMethBytRNAByTissueHeatmap <-
meanMethBytRNAByTissueMat %>%
    Heatmap(
        heatmap_width = unit(5.4,"inches"),
        heatmap_height = unit(11.8,"inches"),
        column_title_gp = gpar(fontsize = 11, fontface = "bold"),
        #"Fetal Tissue tRNA methylation",
        column_title = "Mean Methyation by tRNA\nin Fetal Tissue", 
        name = "Meth",
        column_names_rot = 45,
        row_dend_width = unit(0.2, "npc"),
        row_names_gp = gpar(fontsize = 6),
        row_title = "tRNA gene",
        col = col_fun
    )

png(
    filename = "graphics/meanMethBytRNAByFetalTissueHeatmap.png",
    #width = 9, height = 12, units = "in", res = 192
    width = 6, height = 12, units = "in", res = 192
)
meanMethBytRNAByTissueHeatmap
dev.off()

null device 
          1 

meanMethBytRNAByTissueHeatmap

meanMethBytRNAByTissueHeatmapAAsplit <-
meanMethBytRNAByTissueMat %>%
  Heatmap(
    heatmap_width = unit(5.4, "inches"),
    heatmap_height = unit(11.8, "inches"),
    column_title_gp = gpar(fontsize = 16, fontface = "bold"),
    column_title = "Mean Methyation by tRNA\nin Fetal Tissue",
    name = "Meth",
    column_names_rot = 45,
    row_dend_width = unit(0.2, "npc"),
    row_names_gp = gpar(fontsize = 6),
    row_split = rownames(meanMethBytRNAByTissueMat) %>%
      gsub(pattern = "tRNA-(\\w+)-\\w+-\\w+-\\d+", replacement = "\\1"),
    row_title = "tRNA gene",
    col = col_fun
  )

png(
  filename = "graphics/meanMethBytRNAByFetalTissueHeatmapAAsplit.png",
  #width = 9, height = 12, units = "in", res = 192
  width = 6, height = 12, units = "in", res = 192
  )
meanMethBytRNAByTissueHeatmapAAsplit
dev.off()

null device 
          1 

meanMethBytRNAByTissueHeatmapAAsplit

pseudoSplit <- tRNA$pseudo
names(pseudoSplit) <- tRNA$tRNAname
# pseudoSplit[rownames(meanMethBytRNAByTissueMat)]

meanMethBytRNAByTissueHeatmapPseudosplit <-
meanMethBytRNAByTissueMat %>%
    Heatmap(
        heatmap_width = unit(5.4, "inches"),
        heatmap_height = unit(11.8, "inches"),
        column_title_gp = gpar(fontsize = 16, fontface = "bold"),
        column_title = "Mean Methyation by tRNA\nin Fetal Tissue",
        name = "Meth",
        column_names_rot = 45,
        row_dend_width = unit(0.2, "npc"),
        row_names_gp = gpar(fontsize = 6),
        row_split = pseudoSplit[rownames(meanMethBytRNAByTissueMat)],
        row_title = "tRNA gene",
        col = col_fun
    )

png(
    filename = "graphics/meanMethBytRNAByFetalTissueHeatmapPseudosplit.png",
    #width = 9, height = 12, units = "in", res = 192
    width = 6, height = 12, units = "in", res = 192
)
meanMethBytRNAByTissueHeatmapPseudosplit
dev.off()

null device 
          1 

meanMethBytRNAByTissueHeatmapPseudosplit

heatmaply(meanMethBytRNAByTissueMat)

5.3.2 Boxplots

# custom labeler: (character vector in and out)
tRNA_labeller <- function(value) {sapply(value,
    function(n){
    tRNAmeanMethByTissue %>%
        ungroup() %>%
        select(tRNAname,tChr,strand) %>%
        distinct() %>%
        filter(tRNAname==as.character(n)) %>%
        unlist() %>%
        (function(d) paste0(d["tRNAname"]," (",d["strand"],") ",d["tChr"])) %>%
        return()
    }
)}

plots <- tRNAmeanMethByTissue %>%
    group_by(aa) %>%
        do(plot = 
            ggplot(.,aes(tissue,beta)) +
                #geom_point() +
                geom_boxplot(aes(fill=tissue)) + 
                guides(fill=FALSE) +
                ylim(0,1) +
                facet_wrap(~tRNAname,labeller = labeller(tRNAname = tRNA_labeller)) + #,scales = "free_y"
                geom_hline(aes(colour=tRNAge,yintercept = 0.2)) +
                scale_colour_manual(values = c("TRUE"="red","FALSE"="black"),drop=FALSE) +
                labs(   title = paste0("tRNA genes which are ",.$aa," isoacceptors"),
                        subtitle = "mean methylation across all probes in a tRNA gene for each fetal tissue"
                ) + 
                theme(axis.text.x = element_text(angle = 90,vjust = 0.5))
        )

plots$plot

[[1]]

[[2]]

[[3]]

[[4]]

[[5]]

[[6]]

[[7]]

[[8]]

[[9]]

[[10]]

[[11]]

[[12]]

[[13]]

[[14]]

[[15]]

[[16]]

[[17]]

[[18]]

[[19]]

[[20]]

#plots$plot[[3]]
#plots %>% filter(aa=="iMet") %>% pull(plot)

nil <- plots %>% do(out=ggsave(plot = .$plot,
                    filename = paste0("graphics/epiTOCfetalTissuestRNAMethylation/tRNAmeanMethByTissue_" , .$aa , ".png"),
                    width = 12,
                    height = 6.75
                ))
rm(nil)

6 Differential methylation by tissue

6.1 Kruskal-Wallis test

nestedBytRNA <- tRNAprobesMethData %>% 
    nest(-tRNAname)

KWtissueBeta <- nestedBytRNA %>% 
    group_by(tRNAname) %>%
    mutate(test=purrr::map(data, ~kruskal.test(tissue ~ beta, data = .)))#aov

KWtissueBetaGlance <- KWtissueBeta %>% unnest(test %>% purrr::map(broom::glance))
KWtissueBetaGlance %>%
    ggplot(aes(tRNAname,p.value)) + 
        geom_col() + 
        geom_hline(yintercept = 0.05,colour="red") +
        theme(axis.text.x = element_text(angle=90,vjust = 0.5))

No tRNAs for which methylation was significantly different between tissues by Kruskal-Wallis rank sum test

7 Codons

7.1 Mean meth by codon

tRNAmeanMethByTissue %>%
    ggplot(.,aes(tissue,beta)) +#,group=aa
        stat_summary(fun.data = "mean_cl_boot")+
        guides(fill=FALSE) +
        ylim(0,1) +
        facet_wrap(~codon) + #,scales = "free_y"#labeller = labeller(tRNAname = tRNA_labeller)
        geom_hline(aes(colour=tRNAge,yintercept = 0.2)) +
        scale_colour_manual(values = c("TRUE"="red","FALSE"="black"),drop=FALSE) +
        # labs( title = "x",#paste0("tRNA genes which are ",.$aa," isoacceptors"),
        #       subtitle = "mean methylation across all probes in a tRNA gene for each fetal tissue"
        # ) + 
        theme(axis.text.x = element_text(angle = 90,vjust = 0.5))

tRNAs with certain codons seem more methylated than others, some of the more methylated codons seemed to correspond to those that divered from expectation based on frequency of occurance in the exom see below. e.g. GTT, CTT. Therefore I took a look the relationship between mean methylation levels on tRNAs with a given codon and residuals for the relationship between number of tRNA genes and codon fequency.

revTrans <- function(x){
        ntPairs <- data.frame("sense" = c("A","T","G","C","N"),
                              "antisense"=c("T","A","C","G","N"),
                              stringsAsFactors = FALSE)
        strsplit(x, NULL) %>%
        lapply(rev) %>%
        lapply(function(x) ntPairs$antisense[match(x,ntPairs$sense)]) %>%
        sapply(paste, collapse="")
    }

aaBycodon <- tRNA %>% select(aa,codon) %>% unique()
aaBycodon$codon <- aaBycodon$codon %>% revTrans()

AAUsageFreqHuman <- inner_join(
CodonUsageFreqHuman,
aaBycodon,
by="codon") 

NtRNAbyaa <- tRNA %>% 
    group_by(aa) %>% 
    ###
    summarise("NtRNA"=n())

NtRNAbyCodon <- 
tRNA %>% 
    group_by(codon) %>% 
    summarise("NtRNA"=n())


codonUsageFreqNtRNA <- 
inner_join(
    AAUsageFreqHuman %>% group_by(codon),
    NtRNAbyCodon %>% select(codon,NtRNA),
    by="codon"
)

codonUsageFreqNtRNAplot <- 
ggplot(codonUsageFreqNtRNA,aes(countPerThousand,NtRNA)) +
    geom_label(aes(label=codon,fill=aa),size=4,colour="white",alpha=0.8) +
    scale_fill_manual(values = palette) +
    geom_smooth(method = "lm",se=FALSE) +
    labs(   x = "Codon Usage / count per thousand",
            y = "Number of tRNA genes"
    ) + 
    guides(fill=guide_legend(title="Amino Acid")) + 
    geom_text(x=35,y=17,
            label=paste0(   'italic(r)~',
                            '"="~',
                            sprintf("%0.2f",
                                round((cor(codonUsageFreqNtRNA$countPerThousand,codonUsageFreqNtRNA$NtRNA)),
                                digits = 2))
                            ),
            parse = TRUE
    )

codonUsageFreqNtRNAplot

#ggsave(codonUsageFreqNtRNAplot,filename = "graphics/codonUsageFreqNtRNAplot.png",width=8,height=4.5)

codonUsageVsNtRNA <- lm(codonUsageFreqNtRNA$countPerThousand~codonUsageFreqNtRNA$NtRNA)
codonUsageVsNtRNA %>% summary()

Call:
lm(formula = codonUsageFreqNtRNA$countPerThousand ~ codonUsageFreqNtRNA$NtRNA)

Residuals:
    Min      1Q  Median      3Q     Max 
-14.362  -4.088  -1.052   3.631  20.084 

Coefficients:
                          Estimate Std. Error t value
(Intercept)                12.9939     1.7180   7.563
codonUsageFreqNtRNA$NtRNA   0.2020     0.1284   1.573
                          Pr(>|t|)    
(Intercept)               7.92e-10 ***
codonUsageFreqNtRNA$NtRNA    0.122    
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 8.045 on 50 degrees of freedom
Multiple R-squared:  0.04714,   Adjusted R-squared:  0.02808 
F-statistic: 2.474 on 1 and 50 DF,  p-value: 0.1221

model <- left_join(
    augment(codonUsageVsNtRNA) %>%
        rename(
            "countPerThousand" = "codonUsageFreqNtRNA.countPerThousand",
            "NtRNA" = "codonUsageFreqNtRNA.NtRNA"
        ),
    codonUsageFreqNtRNA,
    by=c("countPerThousand","NtRNA")
)

meanMethBycodon <- tRNAprobesMethData %>% 
    group_by(codon) %>% 
    summarise(beta=mean(beta)) %>% 
    mutate(codon=revTrans(codon))

modelMeanMeth <- left_join(model,meanMethBycodon,by="codon")
modelMeanMeth

FetalMeanMethVsCodonUsageNtRNAGap <- 
modelMeanMeth %>%
    ggplot(aes(.resid,beta)) + 
        geom_label(aes(label=codon,fill=aa),size=4,colour="white",alpha=0.8) +
        scale_fill_manual(values = palette) +
        geom_smooth(method="lm") +
        labs(   title="Fetal tissue tRNA methylation",
                subtitle = "Relationship of methylation level to 'gap' between tRNA gene number and usage frequency",
                x = "Residuals from: Codon Usage Count ~ Number of tRNA genes",
                y = "Mean methylation across probes in tRNA genes /proportion methylated"
        )+
        guides(fill=guide_legend(title="Amino Acid")) + 
        geom_text(x=10,y=0.6,
                label=paste0(   'italic(R^2)~',
                                '"="~',
                                sprintf("%0.2f",
                                    round(glance(codonUsageVsNtRNA)$r.squared,
                                    digits = 4))
                                ),
                parse = TRUE
        )

# FetalMeanMethVsCodonUsageNtRNAGap %>% 
#   ggsave(filename = "../graphics/FetalMeanMethVsCodonUsageNtRNAGap.png",
#       width = 12,
#       height = 6.75
#   )

FetalMeanMethVsCodonUsageNtRNAGap

8 Session Info

sessionInfo()

R version 3.6.1 (2019-07-05)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 19.10

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.3.7.so

locale:
 [1] LC_CTYPE=en_GB.UTF-8      
 [2] LC_NUMERIC=C              
 [3] LC_TIME=en_GB.UTF-8       
 [4] LC_COLLATE=en_GB.UTF-8    
 [5] LC_MONETARY=en_GB.UTF-8   
 [6] LC_MESSAGES=en_GB.UTF-8   
 [7] LC_PAPER=en_GB.UTF-8      
 [8] LC_NAME=C                 
 [9] LC_ADDRESS=C              
[10] LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_GB.UTF-8
[12] LC_IDENTIFICATION=C       

attached base packages:
[1] grid      stats     graphics  grDevices utils    
[6] datasets  methods   base     

other attached packages:
 [1] ComplexHeatmap_2.1.0 heatmaply_0.16.0    
 [3] viridis_0.5.1        viridisLite_0.3.0   
 [5] plotly_4.9.0         rvest_0.3.4         
 [7] xml2_1.2.2           broom_0.5.2         
 [9] forcats_0.4.0        stringr_1.4.0       
[11] dplyr_0.8.3          purrr_0.3.2         
[13] readr_1.3.1          tidyr_1.0.0         
[15] tibble_2.1.3         ggplot2_3.2.1       
[17] tidyverse_1.2.1      testthat_2.2.1      
[19] devtools_2.2.0       usethis_1.5.1       
[21] reprex_0.3.0        

loaded via a namespace (and not attached):
  [1] colorspace_1.4-1    rjson_0.2.20       
  [3] ellipsis_0.2.0.1    rprojroot_1.3-2    
  [5] htmlTable_1.13.1    circlize_0.4.8     
  [7] base64enc_0.1-3     GlobalOptions_0.1.0
  [9] fs_1.3.1            clue_0.3-57        
 [11] rstudioapi_0.10     remotes_2.1.0      
 [13] DT_0.9              lubridate_1.7.4    
 [15] splines_3.6.1       codetools_0.2-16   
 [17] knitr_1.24          pkgload_1.0.2      
 [19] zeallot_0.1.0       Formula_1.2-3      
 [21] jsonlite_1.6        Cairo_1.5-10       
 [23] packrat_0.5.0       cluster_2.1.0      
 [25] png_0.1-7           shiny_1.3.2        
 [27] compiler_3.6.1      httr_1.4.1         
 [29] backports_1.1.4     Matrix_1.2-17      
 [31] assertthat_0.2.1    lazyeval_0.2.2     
 [33] cli_1.1.0           later_0.8.0        
 [35] acepack_1.4.1       htmltools_0.3.6    
 [37] prettyunits_1.0.2   tools_3.6.1        
 [39] gtable_0.3.0        glue_1.3.1         
 [41] reshape2_1.4.3      Rcpp_1.0.2         
 [43] cellranger_1.1.0    vctrs_0.2.0        
 [45] gdata_2.18.0        nlme_3.1-141       
 [47] iterators_1.0.12    crosstalk_1.0.0    
 [49] xfun_0.9            ps_1.3.0           
 [51] mime_0.7            lifecycle_0.1.0    
 [53] gtools_3.8.1        dendextend_1.12.0  
 [55] MASS_7.3-51.4       scales_1.0.0       
 [57] TSP_1.1-7           hms_0.5.1          
 [59] promises_1.0.1      parallel_3.6.1     
 [61] RColorBrewer_1.1-2  yaml_2.2.0         
 [63] memoise_1.1.0       gridExtra_2.3      
 [65] rpart_4.1-15        latticeExtra_0.6-28
 [67] stringi_1.4.3       gclus_1.3.2        
 [69] desc_1.2.0          foreach_1.4.7      
 [71] checkmate_1.9.4     seriation_1.2-8    
 [73] caTools_1.17.1.2    pkgbuild_1.0.5     
 [75] shape_1.4.4         rlang_0.4.0        
 [77] pkgconfig_2.0.2     bitops_1.0-6       
 [79] lattice_0.20-38     htmlwidgets_1.3    
 [81] labeling_0.3        processx_3.4.1     
 [83] tidyselect_0.2.5    plyr_1.8.4         
 [85] magrittr_1.5        R6_2.4.0           
 [87] Hmisc_4.2-0         gplots_3.0.1.1     
 [89] generics_0.0.2      foreign_0.8-72     
 [91] pillar_1.4.2        haven_2.1.1        
 [93] withr_2.1.2         nnet_7.3-12        
 [95] survival_2.44-1.1   modelr_0.1.5       
 [97] crayon_1.3.4        KernSmooth_2.23-16 
 [99] GetoptLong_0.1.7    readxl_1.3.1       
[101] data.table_1.12.2   callr_3.3.1        
[103] digest_0.6.20       webshot_0.5.1      
[105] xtable_1.8-4        httpuv_1.5.2.9000  
[107] munsell_0.5.0       registry_0.5-1     
[109] sessioninfo_1.1.1  

9 References

Nazor, Kristopher L., Gulsah Altun, Candace Lynch, Ha Tran, Julie V. Harness, Ileana Slavin, Ibon Garitaonandia, et al. 2012. “Recurrent Variations in DNA Methylation in Human Pluripotent Stem Cells and Their Differentiated Derivatives.” Cell Stem Cell 10 (5): 620–34. https://doi.org/10.1016/j.stem.2012.02.013.

Yang, Zhen, Andrew Wong, Diana Kuh, Dirk S. Paul, Vardhman K. Rakyan, R. David Leslie, Shijie C. Zheng, Martin Widschwendter, Stephan Beck, and Andrew E. Teschendorff. 2016. “Correlation of an epigenetic mitotic clock with cancer risk.” Genome Biology 17 (1): 205. https://doi.org/10.1186/s13059-016-1064-3.